#include "plant_1DSOS_RDM_input_1Dsin.h"
// --------------------------------------------------
plant_1DSOS_RDM_input_1Dsin::plant_1DSOS_RDM_input_1Dsin(double idt,int iDim,int iL,int iMode, double iE_D,double iE_S,double iE_N)
{  
  cout << "plant_1DSOS_RDM_input_1Dsin::plant_1DSOS_RDM_input_1Dsin()\n";
  // Options ------------------------------------------------------------------
  Dim = iDim;
  if(Dim == 1){
    Lx = iL;
    Ly = 1;
    N_CL = 3;
  }
  else{
    Lx = iL;
    Ly = iL;
    N_CL = 5;
  }

  mode = iMode;
  assert(2*iMode<=iL);
  E_D = iE_D;
  E_S = iE_S;
  E_N = iE_N;

  // BASIC BLOCK ------------------------------
  // Input#0 = W (layer/s)
  // Input#1 = A (layer/s)
  // Input#2 = Freq
  // Input#3 = T (K)
  int idim[4] = {1,1,1,1};
  // Output#0 = h (layer/s)
  // Output#1 = r2 (layer^2)
  // Output#2 = m2
  // Output#3 = a (site)
  // Output#4~#5 = values of eigenfunction phi_1 ~ phi_2
  int odim[6];
  odim[0] = 1;
  odim[1] = 1;
  odim[2] = 1;
  odim[3] = 1;
  odim[4] = (mode+1); // The mode with sin*sin
  odim[5] = (mode+1);
  block::initialize(4,idim,6,odim,idt);

  lattice = new MCArray(Lx);
  cout << "MCArray created successfully\n";
  //	cout << "block::input_dim[]=" << block::input_dim[0] << ", " << block::input_dim[1] << "\n";
  // States -------------------------------------------------------------------
  // lattice = new int*[Ly];
  //	lattice = new MCArray(Lx);
  lattice_stat =  new int*[Ly];
  xcov_val = new double[Lx];
  for(int i=0;i<Ly;i++){
    // lattice[i] = new int[Lx];
    lattice_stat[i] = new int[Lx];
  }
  // Initialize random number generator -------------
  gsl_rng_env_setup();
  gsl_T = gsl_rng_default;
  gsl_r = gsl_rng_alloc(gsl_T);
  gsl_rng_set(gsl_r,time_seed());

//  double SQRT_TWO = sqrt(2.0);
//  double SQRT_PI = sqrt(pi);
  reset();

  // Save a log -------------------------------------
  ofstream flog("plant_1DSOS_RDM_input_1Dsin.log");
  flog << "Dim     = " << Dim << "\n";
  flog << "Lx      = " << Lx << "\n";
  flog << "Ly      = " << Ly << "\n";
  flog << "E_D     = " << E_D << "\n";
  flog << "E_S     = " << E_S << "\n";
  flog << "E_N     = " << E_N << "\n";
  flog.close();
  cout << "END OF plant_1DSOS_RDM_input_1Dsin::plant_1DSOS_RDM_input_1Dsin()\n";
}

// --------------------------------------------------
plant_1DSOS_RDM_input_1Dsin::~plant_1DSOS_RDM_input_1Dsin()
{
  cout << "Begin of ~plant_1DSOS_RDM_input_1Dsin()\n";
  int i;
  cout << "before delete lattice\n";
  delete lattice;
  cout << "deleted lattice\n";
  for(i=0;i<Ly;i++){
    // delete []lattice[i];
    delete []lattice_stat[i];
  }
  
  //delete []lattice;
  delete []lattice_stat;
  delete []xcov_val;
  delete []p_migr;
  delete []p_desorp;
  gsl_rng_free(gsl_r);
  cout << "End of ~plant_1DSOS_RDM_input_1Dsin()\n";
}
// --------------------------------------------------
void plant_1DSOS_RDM_input_1Dsin::reset(){
  block::time = 0.0;

  // Inputs
  T = 0;
  W = 0;

  // States
  int x,y;
  for(y=0;y<Ly;y++){
    for(x=0;x<Lx;x++){
      lattice[y][x] = 0;
    }
  }
  
  for(x=0;x<Lx;x++){
    //cout << "lattice->SetElement(" << x << ",0.0)\n";
    lattice->SetElement(x,0.0);
  }
  for(y=0;y<Ly;y++){
    for(x=0;x<Lx;x++){
      init_local_stat(x,y);
    }
  }

  // Outputs
  h  = 0.0;
  r2 = 0.0;
  m2 = 0.0;
  a  = Lx;
  cout << "END of plant_1DSOS_RDM_input_1Dsin::reset()\n";
}
// --------------------------------------------------
void plant_1DSOS_RDM_input_1Dsin::update(double sysTime){
  if(time>sysTime)
    return;

  double rdm;
  while(time<=sysTime){
    // Update states
    if(iT!=T)
      set_temperature(iT);
    if(iW!=W)
      set_ads_rate(iW);

    rate_ttl = ads_ttl+migr_ttl+desorp_ttl;
    rdm = RandRand()*rate_ttl;
    if(rdm<ads_ttl)
      adsorption();
    else if(rdm<=migr_ttl+ads_ttl)
      migration();
    else
      desorption();

    time += OneMCTime(rate_ttl);
  }
  
  h2mode();
  outputs[4][0] = getH();
  outputs[5][0] = getR2();
  outputs[6][0] = getM2();

  return;
}
// --------------------------------------------------
double plant_1DSOS_RDM_input_1Dsin::getM2(){

}
// --------------------------------------------------
double plant_1DSOS_RDM_input_1Dsin::getR2(){
}
// --------------------------------------------------
double plant_1DSOS_RDM_input_1Dsin::getH(){
}
// --------------------------------------------------
double plant_1DSOS_RDM_input_1Dsin::getTime(){
}



// --------------------------------------------------
void plant_1DSOS_RDM_input_1Dsin::set_temperature(double T_set){
  T = T_set;
  init_migr_rate();
  init_desorp_rate();
  k_m2d = exp(-(E_D-E_S)*EV2KB/(T+0.0));
  return;
}
// --------------------------------------------------
void plant_1DSOS_RDM_input_1Dsin::set_ads_rate(double W_set){
  W = W_set;
  
  return;
}


// --------------------------------------------------
void plant_1DSOS_RDM_input_1Dsin::add_remove_atom(int x,int num){
}
// --------------------------------------------------
void plant_1DSOS_RDM_input_1Dsin::GetNeighbors(int x,int *xpos){
  
}
// --------------------------------------------------
void plant_1DSOS_RDM_input_1Dsin::adsorption(){
  
}
// --------------------------------------------------
void plant_1DSOS_RDM_input_1Dsin::migration(){
}
// --------------------------------------------------
void plant_1DSOS_RDM_input_1Dsin::migr_site(int *migr_x, int *migr_y){
}
// --------------------------------------------------
void plant_1DSOS_RDM_input_1Dsin::desorp_site(int *desorp_x, int *desorp_y){
}
// --------------------------------------------------
